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ABSTRACT 

A detailed solution of an initial value problem of a vertically localized initial pertur- 
bation in rotating magnetized vertically stratified disk is presented. The appropriate 
linearized MHD equations are solved by employing the WKB approximation and the 
results are verified numerically. The eigenfrequencies as well as eigenf unctions are ex- 
plicitly obtained. It is demonstrated that the initial perturbation remains confined 
within the disk. It is further shown that thin enough disks are stable but as their 
thickness grows increasing number of unstable modes participate in the solution of 
the initial value problem. However it is demonstrated that due to the localization of 
the initial perturbation the growth time of the instability is significantly longer than 
the calculated inverse growth rate of the individual unstable eigenfunctions. 
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1 INTRODUCTION 

Magnetohydrodynamic (MHD) instabilities play a major role in a great variety of astrophysical and spa ce applica t ions. Their 
importance is epitom ized by the magneto-rotational instabili ty (MRI) that ha s been discovered by IVelikhov (1959) and 
Chandrasekhar ( 1960l ) for infinite cylinders and rediscovered by Balbus fc Hawlev (1991) in astrophysical context. It is widely 
believed to be one of the prime candidates to provide a viable clue to solving the age old puzzle of the outwards transfer of 
angular momentum in a plethora of astrophysical disk configurations. 

However, for a case of finite size system, the common practice in stability analysis to expand the state vector of the 
linearized dynamical system (the MHD system of equations in our case) in plane waves as: \i(z,t) = Ae"™ 4 " 1 "'* 2 , where z 
is the spatial coordinate (for simplicity only one spatial coordinate will be considered), k is the wave number, and u is the 
natural frequency of the system, is not applicable especially for inhomogeneous systems. Instead, for such cases the appropriate 
bounda ry value proble m (BVP) for a given set of initial conditions (IC) should be solved (see for example the well-known 
work of Landau (194^)). The solution of such problems is facilitated by obtaining the natural frequencies of the bounded 
system uj n for which the state vector of the linearized dynamical system may be written as u n (z, t) = A n (g)e~ I "" t , wher e 
A„(z) are the e i genfu n ctions of the BVP sub jected to specific boundary conditions (BC) (see also Sano fc Mivama ( 19991 ); 
Coppi fc CoppI l|200lh ; ICoppi fc Kevesl (|2003T l). The main purpose of the current work is to employ such approach in order 
to investigate the dynamical evolution of localized initial perturbations in rotating magnetized disks of finite thickness. 

In reality, small perturbations, especially in a system of finite dimensions, do not have the form of a single monochromatic 
wave, but are rather a superposition of individual waves, i.e., wave packets. Furthermore, the asymptotic behavior at large 
times of the wave packet may significantly differ from that of any one of its individual components. Thus, under certain 
circumstance, even though some of the components of the wave packet may individually grow with time without bound, the 
wave packet as a whole can remain bounded at a given place and even decay to zero as the packet is convected away. In such 
cases the initial perturbations give rise to what is defined as convective instabilities. In other cases, the exponentially growing 
components of the wave packets may indeed cause the perturbation to grow without bound at each place. In such cases, the 
perturbations are defined as absolute instabilities. Furthermore, the BC could significantly change the results of the stability 
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analysis (global stability). It is therefore of utmost importance to investigate the asymptotic (in time) development of the 

solutions of the appropriate Cauchy problems in order to determine whether the system is stable or not. 

V a rious sets of rule s in order to distinguish between convective and abso l ute instabilities have been gi v en bvlLandau fc Lifshitz 
lll959h:ISturrockl lll95gf):lFainberg et all ljl96lh : lAkhiezer fc Polovinl (jl97lh : lLifshitz fc Pitaevskiil |l98lh : IPrazinl [|2002h . while 



Huerre fc Monkewitzl |l990h have reviewed more recent developments of the theory pertaining to hydrodynam ic stabil i ty. Th e 



importance of studying the influence of a finite system size on stability analysis has been demonstrated by iBudkerl (|195 ! 
Sturrock who reported that the size of system could play a stabilizing role for the two stream instability. 



One sh ould note also that such understanding may play an im portant role in recent attempts to observe the MRI in the 
laboratory ijNoguchi et alJbooi iRudiger et allbood ; |ji et alJbood ). 



2 THE MAGNETO ROTATIONAL INSTABILITY 



2.1 Formulation of the problem and formal solution 

The effect of the finite size of the system is investigated by re-examining the MRI for the case of rotating disks with finite 
thickness. The basic MHD equations that describe the dynamical development of the system are: 



p%- = -VP + -J x B + pG, 

at c 



V • V = 0, 



dB_ 

dt 



-cV x £, VxB 



47T 



-J, V-B = 0, 



(1) 



(2) 



where G is acceleration due to gravity, and c is speed of light and the rest of the variables have their usual meanings. The 
expression for electric field E in ideal plasmas is given by: 



E = --V x B. 

c 



(3) 



In order to simplify the calculations (and with no loss of generality of the results) we assume that the Brunt- Vaisala frequency 
is small in comparison to all other characteristic frequencies in the system and can hence be neglected. Thus, the linearized 
MHD system of equations that describe radially independent perturbations in a Keplerian disk under the influence of a 
constant axial magnetic field, may be represented in the following way: 



du 

~dt 



dz 



)u = 0, 



where u(z,t) = (v r ,v tfi ,b r ,b v> )' 1 
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(4) 



(5) 



where p(z) is the steady state density profile, and v r ,v v ,b r , and b v are the perturbed radial and azimuthal velocities 
(v r ,v v <§C f2r), and the radial and azimuthal components of the magnetic field (b r ,b v B z ), respectively. In eq. (0 the 
density is normalized to its value at the disk midplane po — p(z = 0), the velocities are normalized to \f2c s (c s is the sound 
velocity), the time is normalized to the inverse steady state angular velocity Q, the magnetic field is scaled with the steady 
state axial magnetic field B z . Finally, lengths are scaled by t/JIVa/Q where Va = B z /^/4npo is the Alfven velocity and /3 is 
the familiar plasma parameter given by (5 = 2c 2 s /Va- 

As a specific example consider a thin isothermal Keplerian disk. In order to model the finite thickness of the disk the 
following profile 

p{r,z) = p (r) exp(-z 2 ) 

is assumed for the density. Due to the independence of the perturbations on the radial direc tion r is from now on m erely a 
parameter. This is tantamount to the local approximation in r with k r = 0, which according to lBalbus fc Hawlevl 1)19911) is the 
most unstable case. For that case, assuming that the IC are u(z,0) = (v r (z, 0), 0, 0, 0), and employing the Laplace transform, 
the set of eqs.Q is reduced to the following single inhomogeneous ordinary differential equation (ODE): 



L[b v ] = 2iujf3 2 v' r (z,Q) 



(6) 



where 
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A = \fpu) represents the re-normalized spectral parameter w, and b v is the Laplace transform of the azimuthal component 
of the perturbed magnetic field. The solution of eq.© should satisfy the following BC: 

M±oo) = 0, (e z \)J i=±oo = -2iyfp\b r (±oo) - A 2 M±oo) = 0. (8) 

Boundary conditions © are obtained due to the requirement that at z — * ±00, where p — > and B z — const, the energy 
flux of the perturbation is finite. 

The solutions of the BVP ©-(JS]) can be constructed by the countable set of the solutions of the homogeneous equation 
L\b v \ = 0, subject to homogeneous boundary conditions. The solutions of that problem are termed eigenfunctions of the 
homogeneous equation. Assuming completeness and orthogonal property of the eigenfunctions allow to express the solution of 
eq.([6} in a simple way for an arbitrary IC by using Green's function (see appendixes A and B for proof of orthogonality, and 
construction of the Green's function). The latter is represented by an appropriate set of eigenfunction b Vl n{z). Thus, carrying 
out the inverse Laplace transform, the solution of eq.(j4]| for t ^ is given by the following expression: 

M ^)-4£« t . >e( ^ w *,"'<' t '--<".)^<"r<»»„ m 

n 

where b ip . n (z) are the eigenfunctions of the homogeneous equation, d n are the coefficients of the expansion of 1^(2,0) i n 
the complete set b v<n (z), and w±(n) are the natural frequencies that correspond to the eigenfunctions b v , n {z). 



2.2 WKB solution 



We turn now to obtaining asymptotic solutions for b lP]n (z) in the limit A > 1 [see IjErokhin fc MoiseevNl973h as an example 
of employing the WKB approximation for some problems in inhomogeneous plasmas]. In that case the eigenfunctions may be 
represented by the following version of the WKB approximation: 



b v ,k = e~ z Mz) = e"" e Sfc(z) , k = 1, . . . ,4. 

Inserting expression (|10p into the homogenous part of eq.® yields: 

2S" 2 + S" 



S' 4 + A 2 <^ 3P(z, X)S + 2Q(z, XjS' + R(z, A) + 3S" 



P(z,\) + 



where 



A^ 



+ 4S'S"' + S IV 



(10) 



(11) 



e~ z (2 + 3/3/A 2 )- (6 + 4z 2 )/A 5 



and 



R(z,X) = e~ 2z (\ 2 -f3) 
By keeping the leading terms in A, eq.|TT} can be written as 
S' 4 + A 2 [3P(z, \)S' 2 + 2Q(z, X)S' + R(z, A)] = 



Q(z, A) = — z 
2e _z2 (2 + 3/3/A 2 ) 



'(2 + 3/3/A 2 )+4/A 2 



(12) 



The solutions of eq. (|12p provide valid approximations for the eigenfunctions of L throughout the disk except in the vicinity 
of the turning points. Assuming that the latter occur at \zq\ 3> 1, it can be shown that the last term in eq. (|12[) becomes 
exponentially small in comparison with other terms (this assumption will be validated later on). In that case the equation 
that determines the turning points is given by 



H^ = A 2 P 3 (2o,A) + Q 2 (2 ,A) = 0. 



Asymptotic solution of eq. (|13[) for large A yields the following expression for the turning points: 

I 2A 2 
20 = ± V ln In^ +O ^ 3 )' 

as well as the following form of the discriminant in the vicinity of the turning points: 



D(«,A) mK(z- 



ZQ 



K 



dD{z,\). 

Tz ' 



-6z (l + 52 ) 



(13) 



(14) 



(15) 



In the outer regions namely z < —zo and z > 20 all four functions Sk(z) are real, whereas in the inner region —20 < 2 < 20 the 
solutions for Sk(z) acquire non zero imaginary part. Thus, the turning points separate the inner range —20 < 2 < 20 where 
the eigenfunctions are oscillatory from the outer regions 2 < —zq and 2 > 20 where the eigenfunctions decrease exponentially 
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Figure 1. The solutions of eq. l|12| l. obtained numerically for = 10, A = 25. Each of the dashed lines (A, B) represents the same real 
part of two complex conjugate solutions. 



away from zq and have no zeros at a finite distance z. The behavior of the various solutions of eq, (|12p as well as the turning 
points may be seen in Fig.l. 

It is obvious from assumption (|10|) and BC |[SJ that only solution 1 is admissible. The turning points are obtained as 
the merging of solution 1 and 2. In the region between the turning points two solutions are admissible, that are complex 
conjugates of each other. Again, due to BC ([8} those solutions are marked by A that depicts their common real part. 

In the vicinity of the turning points the solutions assumed in eq. (|10p are invalid and instead an appropriately re-scaled 
eq.© is solved. That solution is then matched asymptotically to the outer solution given by (|10[) and (|12[) where S' is obtained 
from branch 1 in Fig.l. 

In order to obtain an asymptotic expression for S" in branch 1 it is noticed again that R(z, A) is asymptotically much 
smaller than the rest of the terms that multiply A in eq. (|12[) . As a result, one of the roots is asymptotically zero (branch 3) 
while the rest three roots are obtained from the Cardano solution of the reduced eq. (|12|l . Thus, the root that is represented 
by branch 1 (the only admissible root) gives rise to the following solution in the outer region close to the turning point: 

C 1 

■ exp [ 



y/K(zo - m\ (16) 



In the inner region the solution near the turning points is given by 
1 



\K\(z - z)- 



Ci exp [i j yjK{£-za)d£} + C 2 exp [-i j V - *o)d£] j . (17) 



The coefficients (Ci,2) are determined by the connection formulas tha t express the asymptotic matching of the solution found 



mp 

to the left ()17p and to the right ()16|) of a turning point zo . Following lLandau fc Lifshita ( 19771 ) the connecting relations are 



Ci, 2 = ^exp(ii^). 

Applying the same rule to the region close to the left turning point located at — zo and requireing that the two expressions 
are the same throughout the region —zo < z < zo (the sum of their phase must be multiple of 7r) results in the following 
expression for the eigenfunction: 



b v , n (z,u±{p,n)) = 



cos 



K<&(z) — — 



(18) 



where &(z) = f*° e~~d(^. The analytical expression (|18p for the solution in the inner region has been obtained by noticing 
that away from the turning points S'Q(z,X) is smaller than the rest of the terms in eq. (|12|l and hence to leading order the 
latter is a biquadratic equation. The discrete set of eigenvalues uj is now determined from the Bohr-Sommerfeld relation that 
is given by: 
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K = J | ^3 + 2w±(/3, n)» + ^9 + 16w±09, nf = ^^-y (n + 5) , (19) 

The integer n in the last expression describes the number of zeros of the solution within the inner region [see eg. (|18f) ] . whereas 
w± {0, n) are the natural frequencies of the system that correspond to such eigenfunctions. Strictly speaking expression (|19|l 
is valid for n >> 1. Nevertheless, as is often the case for WKB-type solutions, numerical calculations indicate that expression 
(|18[) provides a close approximation also for n's as low as 1 (see Fig. 2). It is also evident from Fig. 2 that the wave length of 
the perturbation depends on z and approaches the disk thickness close to the turning points. 

As n — is the first excited unstable mode, it is easy to show that the disk is stable for ^ 0* where 0* = 
^/(24[£r/(^)] 2 ) (n = 1,0. = 3 7 r/(8[ J Br/(^)] 2 )). 



2.3 Numerical solution 

To further follow the development of a localized initial perturbation governed by eq. (|4]), consider for example an IC given as: 

v r (z,0) = e~ z2/A \ (20) 

and zero for the rest of the initial values of the perturbed physical variables. It is easy to see that the perturbation is localized 
within the disk, if A << 1. The dynamical development of the localized initial perturbation that is described by eq.@ depends 
on the disk "thickness" defined by the plasma parameter (notice that now velocities are scaled by the Alfven velocity which 
means that the disk thickness is 2y/0~ ). If the disk is thin enough such that 0^0* all the natural frequencies w±(n) that are 
obtained from eq. (O are real and hence the disk is stable. In that case following <(9j the profile of the perturbed azimuthal 
magnetic field represented by two identical wave packets that move in opposite directions namely up and down. Upon reaching 
the upper and lower turning points (located at ±Zo) the two packets are reflected and continue their motion through the disk. 
The amplitudes of the identical wave packets change in time like oc sin(f). It should be noted however that this picture is 
valid if the initial perturbation is localized enough. As example for such behavior, b v (z,t), which has been obtained by using 
eq.([9} A = 0.1, and = 3ir/8 is depicted in Fig. [3] In this particular case all eigenfunctions have real natural frequencies 
[w±(l) = 1.54, 3.64, u>± (3) = 3.99,6.01,...] and there are no growing modes at all (the wave packets indeed do not grow in 
time, just move back and forth inside the disk with local Alfven velocity). 

For "thick" disk (0 > 0*) some of the natural frequencies have an imaginary part. As the thickness of the disk increases 
so does the number of unstable modes. Thus for long time b lp (z,t) ~ b v , tn (z)exp(\ Qui_(n) | t) where n corresponds to natural 
frequency with a biggest imaginary part. However, for sufficiently localized initial perturbation the fastest growing mode has 
sufficiently small initial amplitude. The perturbation then keeps it's original shape for a long time until a growing mode 
emerges out of the wide initial spectrum and is of order of the amplitude of initial wave packet. Raising the value of to 10, 
one unstable mode, with growth rate 7 = 0.749 enters the spectrum of the eigenmodes. Employing again eq.© together with 
eigenfunctions b V: n{z) results in Fig.4. It is indeed seen that after long time the unstable mode dominates the perturbation. 
Thus, the latter eventhough starting as a well confined perturbation develops into a global instability. It should be noted 
however, that the time it takes the unstable mode to emerge from the wide spectrum that makes up the localized initial 
perturbation is about 4 times longer than the predicted growth time (inverse growth rate). Indeed, as the initial perturbation 
is more localized so the growth time of the most unstable mode is increased relative to its linear predicted value. 
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Figure 3. Time evolution of the profile of the perturbed azimuthal component of the magnetic field. The different profiles were calculated 
at times [0.3, 1.0, 1.5] for (3 = 3tt/8, A = 0.1 . 




Figure 4. Long time evolution of the profile of the perturbed azimuthal component of the magnetic field for /3 = 10, A = 0.05. The 
different profiles were calculated at times [0.3, 1.0, 6.0]. 



It is finally instructive to plot the wave number n* of the most unstable mode as a function of /3 (Fig. 5). Asymptotic 
estimation reveals that n* is proportional to yf]3. This implies that the wavelength of the most unstable mode is of the order 
of h/n* ~ Va/Q (where h is the thickness of the disk) which to leading order does not depend on f3. 
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Figure 5. Axial standing wave number of the most unstable mode as a function of the plasma beta. 
2.4 Limit of infinite "thickness" 

It is instructive to examine the limit of thick disk such that effect of boundaries may be ignored. Such an assumption is 
justified when the Alfven traverse time is longer than the time of interest (say inverse growth rate). Renormalizing therefore 
lengths to Va/& solution is sought for z <C y 7 /?, where 2y/~j3 is now the disk's thickness. It is readily seen that for such case 
btp, n (z) —* b v (k,z) = cos fez for odd initial perturbations [see ea. (|18p where k is defined as yj ir /2/3(n + 1/2)] so that the 
eigenfunctions can be approximated by plane waves. Thus, due to l/yffi — * 0, the sum ([9]) is replaced by an integral over k 
and eq.© can be written as 



A{k)e lkz dk 
D(k,u) 



(21) 



where C is the contour of integration in the complex u) plane that is a straight line parallel to the real axis and passing above 
all singular points of the integrand, A(k) is the Fourier transform of the right hand side of eq. ((6} and the denominator of 
the integrand is given by 



D(ui, k) = [oj — k ) — lu 



3fc 



(22) 



It is important to notice now that since the initial perturbation is well localized in space A(fe) in contrast is close to a constant 
that is proportional to the localization length. As a result it can be shown now tha t b^(z,t) obtained by using (|21[> is not zero 
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only if the integrand has a singularities of order two (branch points) (see in details in lFainberg et al 
( 19711 )). Thus, examining the equation D(uj, k) = [the classical MRI dispersion equation IjBalbus fc Hawlevl ll991)] it is clear 



that there is just one branch point that gives rise to instability, which is uo = 3i/4, k = 
behavior of the perturbation is given by 

b v {t) ~ e 3t/4 



15/4. Consequently, the long time 



(23) 



As remarked above this result is due to the fact that the initial perturbation is localized in space. If, however, (a not too 
physical) monochromatic perturbation is considered, A(k) is a delta function and the familiar purely exponential growth is 
recovered. 



3 CONCLUSIONS 



The imp ortance of solvin g the initial value probl em wi th some appropriate boun dary conditions is highlighted. In the classical 
works of lVelikhovl |l959h : Ichandrasekhari l|l960l ). and lBalbus fc Hawlevl |l99lh . an infinite cylinder has been considered and 
the effects of the boundary conditions were neglected. Therefore, in those works, naturally, the thickness of the disk does not 
play any role and consequently cannot influence the extent of the domains of instability. In spite of that, physical intuition and 
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insight have led lBalbus fc Hawlevl (| 199 ll ) to conclude that j3 w 1 is the lower limit for the disk thickness for the occurrence of 
the classical MRI in disks as for smaller values of /3 the thickness of the disk is smaller than the wave length of the dominant 
unstable mode. In the current work the stabilizing effects of the boundaries are taken explicitly into account and hence the 
threshold thickness of the disk may be calculated easily from the WKB solutions obtained in Sec. 2. Furthermore, the number 
of unstable modes as a function of f3 may be estimated with the aid of that solution. Thus, for example, it may be shown 
that there are only three unstable modes within the disk as long as the its beta value is less than ~ 15. Such information 
is significant for the study of the consequent nonlinear development of the instability. In addition, the shape of the unstable 
(as well as the stable) perturbations has been obtained explicitly and as may be seen by its expression (|18|l may differ 
significantly from the plain waves assumed in the classical works on MRI's. It is finally interesting to note, that going to the 
limit of thick disks (see Section 2.4) and assuming a localized perturbatio n results in a reduced r ate o f growth in comparison 
to the less realistic monochromatic classical result [see Eq.(23J|]. Recently. Icoppi fc Coppl l|200ll ), and lCoppi fc Kevesf (|2003h 
have considered the effects of the axial localization of the perturbations by studying axisymmetric ballooning modes in finite 
disks. Such modes are characterized by finite values of the radial wave vector k r . The appropriate turning points within the 
disk were found and consequently a discrete spectrum of eigenfunctions and eigenvalues was obtained. The growth rates were 
found to be smaller than their "long-cylinder" counterparts. In particular, it is shown there that the k r — > limit cannot be 
obtained from their scaling s and asymptotics, as no turning points exist within the disk in that limit. In that sense the present 
work is complementary to (jCoppi fc Kevesll2003h as it analyzes that very limit k r — > 0. Indeed, under the current scaling and 
asymptotic expansion turning points are found within the disk and the corresponding discrete spectrum is obtained. 

The linear MHD equations have been employed in order to study the stability and time behavior of a rotating stratified 
Keplerian disks whose density decreases with height. The full solution of the dynamical evolution in time of an localized 
Gaussian wave packet is explicitly derived and it's long ti me behavior is discusse d. It is proven analytically (WKB) that MRI 
can be suppressed in sufficiently thin disk /3 < /?» ~ n /ll |Miller fc Stonel (|2000l ) observed from numerical simulations that if 
P « 1 the disk becomes MRI stable). Nume rical solutions of the Ca uchy problem however, indicates that /3* is of order 1. This 
result is consistent with values obtained bv lSano fc Mivama (1999). For thicker disks, the number of discrete unstable mode 
increases with the thickness of the disk. However, due to the localization in real space the initial amplitude of the unstable 
mode is diminishingly small and the time it outgrows the original wave packet is significantly longer than its predicted inverse 
growth rate. Thus, it has been demonstrated numerically the growth time of a perturbation that is initially localized within 
the inner 5% of the disk, may be an order of magnitude longer than the inverse growth rate of the fastest growing unstable 
mode. 

The considering of the boundary effects may have a large impact on the design of laboratory experiments to model MRI, 
where a magnetic field have to be quite strong and devices cannot be very high. 
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APPENDIX A: ORTHOGONALITY 



Both sides of the homogenous equation L[b ip , n (z)] = are first multiplied by b Vim {z) and then integrated over the interval 
(— oo,oo). Integrating by parts and setting the boundary terms to zero due the homogeneous BC © yields 



b<p,n 

) 

Defining 



(e z (e z b'^m) } + 3/3 (e z b' Vt ^\ 



dz + 2A„ 



dz + A 2 ,(A 2 - P) I b v , m b v , n dz = (Al) 



btp,rnbip,ndz, Q'2,, 



dz = 



( z2 7? V 



dz. 



and noticing that to each eigenfunction b v , n correspond two eigenvalues A^ x 2 allows rewriting (|A1[I in the form of the following 
set of homogenous linear equations: 

Xn,l(^n,l — P) — Am(A^ — 0) 2(A^ ljl — A^ ) 
-^n,2(An,2 — P) — A m (A m — 0) 2(A Uj2 — A m ) 



AQ — 0, Q = (Ql,m,ri, Q2,m,n) 



(A2) 



Implying that different eigenfunctions (n 7^ m) have different eigenvalues results in: 
det(A) = 2(A^ - A^ 2 )(A^! - A^A 2 ,, - \ 2 m ) ± 0. 
This implies that 

Ql,m,n — Q2,m,n — 0. 



(A3) 



(A4) 



In the opposite case of identical eigenfunctions (m = n) the equality det(A) = means that Qi. n , n 7^ 0. Combining the latter 
with eq. ([A4[) gives 



(6, 



(A5) 



APPENDIX B: GREEN'S FUNCTION 

To obtain the expansion coefficients c n of a solution b v {z) for the inhomogeneous equation L[b^(z)] — —f(z) in terms of the 
eigenfunctions b^^iz) the solution expansion b lp (z) = cib v ,i(z) and expansion of the r.h.s of the inhomogeneous equation 
/(z) = d m b Vim (z) is substituted into the inhomogeneous equation. The result is multiplied by b^^iz) and integrated 
over the interval (— 00, 00). The result is 

/oo 
b v , n L[b vJ ] — -d n \\b V!n \\ 2 . (Bl) 
-00 

Defining in addition to Qi, n ,i and Q2,n,i 

"/ — CO 

and employing the homogenous equation = multiplied by b Lp ^ l (z) and integrated over the interval (—00,00) yields: 

Oi.n.i = \K,i\\ 2 S n>l , Q 2 , n ,i = -J^^(A?, 1 + A?, 2 -/3)5„, !) Q 3 ,n,i = ^f- [3/3(A 2 , x - /?) + Af, 2 (2A?, 1 + 3/3)] * n ,,. (B2) 
Thus eq. (|Bl|) can be rewritten as 

Cn(A^i - A 2 )(A 2 , 2 - A 2 ) = -d n , (B3) 
and a solution of the inhomogeneous equation with r.h.s —f(z) has a form: 
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Finally noticing that d n = f(z)b (Ptn (z)dz, a solution of inhomogeneous equation (|B4[) can be rewritten 

II^V, Till °° 



J —00 

where 

G( c) = y K, n (*)KAO g 

„ (A 2 nil -A2)(A^ 2 -A2)||6^ n p 

is a Green's function written in terms of an eigenfunctions b^, n (z). 

This paper has been typeset from a TjtX/ J^Tj^X file prepared by the author. 



